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ABSTRACT 

We utilize ACDM halo occupation models of galaxy clustering to investigate the evolving stellar mass 
dependent clustering of galaxies in the PRIsrn MUlti-object Survey (PRIMUS) and DEEP2 Redshift 
Survey over the past eight billion years of cosmic time, between 0.2 < z < 1.2. These clustering 
measurements provide new constraints on the connections between dark matter halo properties and 
galaxy properties in the context of the evolving large-scale structure of the universe. Using both an 
analytic model and a set of mock galaxy catalogs, we find a strong correlation between central galaxy 
stellar mass and dark matter halo mass over the range Mh a io ~ 10 n -10 13 h~ 1 M Q , approximately 
consistent with previous observations and theoretical predictions. However, the stellar-to-halo mass 
relation (SHMR) and the mass scale where star formation efficiency reaches a maximum appear to 
evolve more strongly than predicted by other models, including models based primarily on abundance¬ 
matching constraints. We find that the fraction of satellite galaxies in haloes of a given mass decreases 
significantly from z ~ 0.5 to z ~ 0.9, partly due to the fact that haloes at fixed mass are rarer at 
higher redshift and have lower abundances. We also find that the Mi/M m ; n ratio, a model parameter 
that quantifies the critical mass above which haloes host at least one satellite, decreases from « 20 at 
z ~ 0 to « 13 at z ~ 0.9. Considering the evolution of the subhalo mass function vis-a-vis satellite 
abundances, this trend has implications for relations between satellite galaxies and halo substructures 
and for intracluster mass, which we argue has grown due to stripped and disrupted satellites between 
0 ~ 0.9 and z ~ 0.5. 

Subject headings: cosmology: theory - cosmology: observations - cosmology: dark matter - galaxies: 

distances and redsliifts - galaxies: clustering - galaxies: halos - galaxies: evolution - 
galaxies: high-redshift - large-scale structure of the universe - methods: statistical 
- methods: analytical 


1. INTRODUCTION 

In the ACDM cosmology, structures form hierarchi¬ 
cally, such that smaller haloes merge to form larger and 
more massive haloes. All galaxies are thought to form 
as a result of gas cooling at the center of the potential 
well of dark matter haloes. When a halo and its ‘central’ 
galaxy are accreted by a larger halo, it becomes a subhalo 
and its galaxy becomes a ‘satellite’ galaxy. In addition 
to mergers, haloes also grow by smooth accretion and 
galaxies grow by in situ star formation when fuel (i.e., 
cool gas) is available. In this paradigm of hierarchical 
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structure formation, there is a correlation between halo 
formation and abundances and the surrounding large- 
scale structure such that more massive haloes tend to 
reside in denser regions (Mo & White 1996). Galaxy for¬ 
mation models assume that the properties of a galaxy are 
determined entirely by the mass and formation history 
of the dark matter halo within which it formed. Thus, 
the correlation between halo properties and environment 
induces a correlation between galaxy properties and en¬ 
vironment. 

The halo model is a useful framework for discussing 
how galaxy clustering depends on the properties of the 
galaxies’ parent dark matter haloes, and it is a useful 
guide for studying the connections between galaxy for¬ 
mation and halo assembly (see Cooray & Sheth 2002; 
Mo, van den Bosch & White 2010 for a review). Halo 
models of galaxy abundances and clustering generally 
consist of the following three types: halo occupation dis¬ 
tributions (HODs; e.g., Seljak 2000; Scoccimarro et al. 
2001; Berlind & Weinberg 2002; Kravtsov et al. 2004), 
conditional luminosity functions (CLFs; e.g., Yang et al. 
2003; Cooray 2006), and (sub)halo abundance match¬ 
ing (SHAMs; Vale & Ostriker 2006; Conroy, Weclisler & 
Kravtsov 2006). Such models are useful for exploring the 
relations between galaxy formation and dark matter halo 
assembly in the context of the large-scale structure of the 
Universe. The stellar to halo mass relation (SHMR) is 
commonly studied one in the literature (Mandelbaum et 
al. 2006; Moster et al. 2010; Behroozi et al. 2010), and 
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one of the goals of this paper is to analyze its evolution FI 
The SHMR and its variants quantify the fundamental 
correlation between the size of central galaxies and the 
parent haloes that host them, spanning from low-mass 
dwarfs and Milky Way galaxies to galaxies in massive 
groups and clusters, which are hosted by massive haloes. 
The ratio has also been used to define the halo mass scale 
of maximum star formation efficiency in galaxies, as it 
provides a ratio of baryons to dark matter as a function 
of halo mass. 

In addition to these relations, we are also interested in 
exploring the distributions of galaxy and halo properties 
with more complex models (e.g., Skibba & Sheth 2009; 
Hearin & Watson 2013; Cohn & White 2013). Moreover, 
sophisticated self-consistent evolutionary models using 
star formation histories and merger trees have been re¬ 
cently developed (e.g., Behroozi et al. 2013c; Yang et 
al. 2013), and these are complementary to semi-analytic 
models as well (e.g., Kang et al. 2012; Q. Guo et al. 2013; 
Contreras et al. 2013; Campbell et al. 2014). 

The literature contains an impressive array of work ap¬ 
plying halo models to measurements of evolving galaxy 
clustering and lensing, and our work is complementary to 
them. For example, galaxy clustering has been modeled 
at a wide range of redsliifts beyond 2 ~ 0.2 in the SDSS 
and BOSS surveys (Wake et al. 2008; White et al. 2011; 
H. Guo et al. 2014); DEEP2 (Conroy et al. 2006; Zheng 
et al. 2007; Watson & Conroy 2013); VVDS (Abbas et 
al. 2010); COSMOS (Leautliaud et al. 2012; Tinker et al. 
2013; McCracken et al. 2015); VIPERS (de la Torre et 
al. 2013); CFHTLenS (Coupon et al. 2015); and Spitzer 
SPT (Martinez-Manso et al. 20150- Other authors have 
analyzed group catalogs and density field reconstruction 
as well. However, previous work at z ~ 1 has been often 
limited to relatively small samples in volumes subject to 
substantial‘cosmic variance’ errors and/or lacked accu¬ 
rate stellar masses or spectroscopic redshifts. 

In this paper, we utilize data from the PRIsm MUlti- 
object Survey (PRIMUS; Coil et al. 2011; Cool et al. 
2013), using volume-limited catalogs constructed from 
a parent sample of over 130,000 galaxies with robust 
redshifts in seven independent fields covering 9 deg 2 of 
sky. In Skibba et al. (2014; hereafter S14), we modeled 
the optical luminosity and color of galaxy clustering at 
0.2 < z < 1.0 with a simple halo model in which the 
HOD parameters are assumed to be constant with red- 
shift. Using new clustering measurements as a function 
of stellar mass and specific star formation rate (Mendez 
et al., in prep.; hereafter M15), we now analyze the scale- 
dependent clustering evolution with improved halo mod¬ 
els to study the evolving relations between stellar mass 
and dark matter halo mass of central and satellite galax¬ 
ies. In order to perform a complete analysis and to obtain 

11 Note that these are sometimes mistakenly referred to as 
‘abundance matching relations’ though they can be inferred with 
a variety of methods: for example, early HOD models produced 
luminosity-halo mass relations (e.g., Peacock & Smith 2000), and 
one could argue that early constraints on galaxy-halo relations 
were obtained with the morphology-density relation (Dressier 1980; 
Postman & Geller 1984) and galaxy kinematics (Zaritsky et al. 
1993; Carlberg et al. 1996). 

12 Details about these surveys can be found at the following 
references, respectively: York et al. (2000); Dawson et al. (2013); 
Davis et al. (2003); Le Fevre et al. (2005); Scoville et al. (2007); 
Guzzo et al. (2014); Heymans et al. (2012); Ashby et al. (2013). 


robust model parameters, we analyze the spatial distri¬ 
butions and abundances of PRIMUS galaxies with two 
types of models: an analytic model and a set of mock 
galaxy catalogs. The former is based on dark matter halo 
statistics and quantities measured from numerical simu¬ 
lations, including the mass function, bias function, and 
density profile. The latter is directly tied to a dark mat¬ 
ter cosmological simulation and halo-finding algorithm. 
These complementary methods are both widely used in 
the literature. 

This paper is organized as follows. In Section [2j we 
briefly describe the PRIMUS galaxy clustering measure¬ 
ments of M15. We describe the analytic halo model of 
galaxy clustering and mock galaxy catalogs in Sec. 0 
and additional details are described in the appendices. 
We then present our HOD model results in Sec. O and 
we discuss the results and their implications, such as for 
the stellar mass-halo mass relation and satellite abun¬ 
dances, in Sec. [5] Finally, we end with a summary of our 
conclusions in Sec. [Gj 

Throughout this paper, we adopt a flat ACDM cos¬ 
mology with H m = 0.27, Ha = 1 — H m , as = 0.8(0, 
n s = 1, and we express units that depend on the Hubble 
constant in terms of h = Hq/100 kms -1 Mpc^ 1 . For the 
halo masses and radii, we define them with a spherical 
overdensity 200 times the mean density of the Universe 
(see Appendix [All . except when specified otherwise. M 
or Mh refers to halo mass and m refers to subhalo mass 
(in units of /i -1 Mq), and Af* refers to stellar mass (M 0 
units). For the stellar masses and star formation rates 
(SFRs), we assume a universal Chabrier (2003) initial 
mass function (IMF). 

2. GALAXY CLUSTERING MEASUREMENTS 

In S14, we presented spatial clustering measurements 
of galaxies with high-quality redshifts in the PRIMUS 
and DEEP2 surveys as a function of luminosity and color 
over the redshift range 0.2 < z < 1.2. In M15, we present 
clustering measurements as a function of stellar mass and 
specific star formation rate (sSFR) over the same red- 
shift range, and these are the measurements we use here. 
The clustering measurements build on the stellar mass- 
dependent clustering of active galactic nuclei (AGN) in 
Mendez et al. (2015a). We quantify galaxy cluster¬ 
ing with two-point projected auto-correlation functions, 
Wp(r p ), which are measured by integrating £(r p ,n) out 
to line-of-sight separations of 7r max = 80 /i -1 Mpc for 
PRIMUS fields and 20 fr 1 Mpc for DEEP2. The cata¬ 
logs are roughly volume-limited, and incompleteness and 
redshift success weights are applied similarly as in S14. 
The correlation function errors are estimated with jack¬ 
knife resampling methods, which are described in S14 
and M15. 

The main properties of the stellar mass-limited sam¬ 
ples are listed in Table [U and the stellar mass and 
redshift limits and distributions are shown in Figure CD 
M15 includes the PRIMUS science fields (CDFS-SWIRE, 
XMM-LSS, COSMOS, DLS, and ESI) as well as the Ex¬ 
tended Groth Strip (Lin et al. 2007) and DEEP2 fields. 
We refer the reader to M15 for details. 

13 Note that these values of and erg are slightly lower than the 
latest cosmological constraints (Planck collaboration et al. 2014). 
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TABLE 1 

Properties of the Stellar-Mass Limited Samples 


sample 

-^min 

<*> 

^max 

logM*, min 

(logM*) 

logM*, max 

log n sa i 

fegal 

Ml 

0.20 

0.49 

0.70 

10.00 

10.65 

12.00 

-2.37i0.06 

1.35i0.02 

M2 

0.20 

0.50 

0.70 

10.50 

10.87 

12.00 

-2.69i0.07 

1.52i0.06 

M3 

0.20 

0.51 

0.70 

11.00 

11.20 

12.00 

-3.45i0.ll 

1.74i0.15 

M4 

0.70 

0.87 

1.20 

10.50 

10.90 

12.00 

-2.7i0.1 

2.13i0.16 

M5 

0.70 

0.91 

1.20 

11.00 

11.20 

12.00 

-3.35i0.07 

2.5i0.3 


Stellar mass threshold samples: redshift limits, stellar mass limits (Mg), and number densities (in units of 10~ 3 h 3 Mpc —3 ), and 
bias. Mean redshifts and stellar masses of the samples are listed, and their median values are not significantly different. The lower 
error bars of the high-redshift biases indicate the results without the COSMOS field (see Appendix ICl and M15 for details). All of 
the values quoted are obtained from M15, except for the number densities and their errors, which are obtained from Moustakas et 
al. (2013; Table 5) and are used here as additional constraints. 



Redshift 


model for constructing mock catalogs. 

Throughout this section, M refers to the mass of dark 
matter haloed (defined using a virial overdensity 200 
times the mean density of the Universe, unless stated 
otherwise) and M* refers to the stellar mass of galaxies. 

3.1. Halo Occupation Distribution 

For these calculations, we use the mean redshifts of the 
observed galaxy samples and their number densities in 
addition to the projected clustering measurements. For 
the analytic model, we assume a Tinker et al. (2008b) 
halo mass function, Tinker et al. (2010) halo bias, and 
Eisenstein & Hu (1998) matter power spectrum. For 
both the analytic model and mock galaxy catalogs, we 
assume a Navarro, Frenk & White (1997) density pro¬ 
file. The details of these quantities are described in Ap¬ 
pendix Q3J 

We use a model of the halo occupation distribution 
(HOD) similar to that in Skibba & Sheth (2009) and 
Zheng et al. (2007), and we refer the reader to these pa¬ 
pers for details. The mean central galaxy HOD is mod¬ 
eled with the following parametrization: 


Fig. 1. - Stellar mass and redshift limits of the five PRIMUS 
samples described in Table [T] Blue and red points indicate star¬ 
forming and quiescent galaxies based on the SFR-M* demarcation 
of Moustakas et al. (2013). 

The stellar masses and SFRs are taken from Mous¬ 
takas et al. (2013), and they are estimated from multi¬ 
wavelength imaging from UV to mid-IR wavelengths us¬ 
ing a Bayesian spectral energy distribution (SED) mod¬ 
eling code (iSEDf it). The fiducial parameters are based 
on the flexible stellar population synthesis models of 
Conroy & Gunn (2010), and exponentially declining star 
formation histories (SFHs) with Gaussian bursts of star 
formation superposed are assumed. Metallicities near 
the solar value (Asplund et al. 2009) are assumed, and 
the time-dependent dust attenuation curve of Chariot & 
Fall (2000) is adopted. The number densities and errors 
in Table [T] are computed from the stellar mass functions 
(SMFs) of Moustakas et al. (2013), which includes a thor¬ 
ough analysis of relevant systematic uncertainties. 

3. MODEL DESCRIPTION 

In this section, we describe the components of our halo 
model of galaxy clustering. We will describe the halo oc¬ 
cupation distribution (HOD) model, the galaxy cluster¬ 
ing model, and the method for constructing mock cata¬ 
logs. Except when stated otherwise, we utilize the same 
model components in both the analytic model and the 


{N m \M,M.) = 


1 + Jm® 

\ 0"logM 


(i) 


where erf is the error function, which assumes a lognor¬ 
mal distribution for the central galaxy conditional SMF, 


^cen (M*| M)dM. 


y/2n ln(10) a\ ogM 


exp 


_ ( log(M*/M«, cen ) 

\ logM 


( 2 ) 

dM* 

~mT' 


The mean satellite galaxy HOD is parametrized by 

( M — A/q \ 

<JV„t|M,M.)= ( ° j , (3) 

and we discuss the satellite conditional SMF in Sec¬ 
tion [572] 

Note that unlike the simple model in S14, in this pa¬ 
per we allow the HOD parameters to evolve with redshift. 
The halo mass (threshold) and stellar mass (threshold) 
are directly related as the M*-Mh a i 0 relation, eri og M 


14 Note that the same formalism may be utilized to model galaxy 
clustering as a function of halo circular velocity (e.g., Umax), rather 
than mass. 























4 


Skibba et al. 


determines the lognormal scatter in the relation, the 
/r = parameter determines the critical mass 

above which haloes typically host at least one satellite 
within the selection limits, and a is the power-law in¬ 
dex of the mass dependence of the efficiency of satellite 
galaxy formation. Mq describes the smooth drop-off of 
the satellite HOD at low halo mass (relative to the mass 
threshold), but in practice a model with ( M/Mi) a pro¬ 
duces similar results (Zheng et al. 20070. 

(N\M) is the mean HOD of all galaxies (i.e., centrals 
and satellites), 

(N ga i\M, M*) = (-/V cen |M, M*)(l + (-V sat |M, M*)) (4) 

The mean galaxy number density, h ga i, is defined in the 
analytic model with the following: 


^ gal ( M% ) 


r M — (M *) dM Mm, z) 

J M min (M*) dM 

x{N cen \M) [l + (N sat \M) 


(5) 


For the central galaxies, we are implicitly assuming 
that their stellar masses have a lognormal distribution at 
fixed halo mass, though it is not known how accurate this 
assumption is (More, Diemer & Kravtsov 2015). For the 
satellite HOD, we are assuming that they follow a Pois¬ 
son distribution (Kravtsov et al. 2004; Zheng et al. 2005). 
We perform most of our model analysis without allowing 
Mo to be a free parameter, as it is poorly constrained and 
we are interested in exploring constraints on more impor¬ 
tant satellite HOD parameters such as Mi = 

The example correlation functions in Figure [2] show 
how the clustering varies with Mi/M m - ln (Eq. [3| over 
the maximum range of the parameter. Note that lower 
values of the parameter imply more satellites and larger 
number densities, but h ga i is primarily determined by the 
halo mass threshold. There is some degeneracy between 
Mi/M m i n and a, because the larger values of the latter 
also increases the clustering strength, though in a differ¬ 
ent scale-dependent manner. In general, when there are 
more satellites hosted by halos of some mass, there will 
be more central-satellite and satellite-satellite pairs and 
therefore stronger two-point clustering. 



Fig. 2.— An example of the dependence of model correlation 
functions on the parameter for a given M m [ n (i.e., as¬ 

suming the same M*-M relation). Clustering predictions of mock 
galaxy catalogs with = 8, 16, and 32 are shown, span¬ 

ning the range of allowed parameter values. Lower values of this pa¬ 
rameter imply a higher satellite fraction and more central-satellite 
and satellite-satellite pairs. 


combined with the occupation distribution of galaxies in 
the halos: 






dM 


dn(M , z) 
dM 


&halo(M, z ) 


(N gal \M) 

^Igal 


( 6 ) 

(Cooray & Sheth 2002; Yang et al. 2003), where the inte¬ 
gration limits are related to the (e.g. stellar mass depen¬ 
dent) selection of the galaxies themselves. In terms of the 
large-scale galaxy correlation function, which depends on 
galaxy stellar mass M*, galaxy bias can be described as 


3.2. Galaxy Clustering Model 

In this section, we describe in more detail how galaxy 
clustering is modeled in the analytic model. 

Galaxies and haloes are biased tracers of the under¬ 
lying distribution of dark matter (e.g., Weinberg et al. 
2004; Zehavi et al. 2011; S14). In the ACDM the¬ 
ory of hierarchical structure formation, the large-scale 
clustering of halos with respect to matter can be de¬ 
scribed with the halo bias parameter: ^hh{r,M,z) ~ 
[frhaio(M, z)f( mm (r, z) (Mo & White 1996; Sheth & Lem- 
son 1999), where the matter correlation function is ob¬ 
tained from the linear or nonlinear power spectrum (Ef- 
stathiou, Bond & White 1992; Smith et al. 2003). In 
the halo model of galaxy clustering, galaxy bias 6 ga j can 
then be inferred from the abundance and bias of halos, 

15 (M/M\ ) a is the way the mean satellite HOD has been tradi¬ 
tionally modeled, so we write M 1 in © because of the Mq param¬ 
eter. 


£ss( r |M*, z) = [6gal(M* 5 Z)\ Cmm (r,z). (7) 

In this paper, we perform a more detailed analysis us¬ 
ing the full spatial correlation function of galaxies at pro¬ 
jected separations of 0.1 < r p < 30 h -1 Mpc (i.e., we are 
not limiting the analysis to large scales as was done in 
S14). Following Skibba & Sheth (2009), we perform our 
halo model calculations in Fourier space. The two-point 
correlation function is the Fourier transform of the power 
spectrum 

dk k 3 P(k) sin kr , N 

T^f 2 — (8) 

In the halo model, P(k) is written as the sum of two 
terms: one that arises from galaxies within the same halo 
and dominates on small scales (the 1-halo term), and the 
other from galaxies in different halos which dominates on 
larger scales (the 2-halo term). That is, 

P(k) = P lh (k) + P 2h (k), 


(9) 
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where, 


Pih{k\M.) 


P2h(k\M *) 


/■M max (M„) dn(M) 

/ dM—LJ(7V cen |M) 

2{N sat \M) u gal (k\M) 

n gal 

, Wat(iV S at-l)|M) Ugal (fc|M) 2 

' - 9 

n gal 

/■M max (M.) dn(M) 

/ dM {N cen \M) 

'M min (M,) 


dM 

1 + (AU|M) %a i(fc|M) 

TT-gal 


b(M) 


,( 10 ) 

( 11 ) 


-Plin(fc), 


and u ga i(/c|M) is the Fourier transform of the galaxy den¬ 
sity profile, which is closely related to the subhalo and 
dark matter density profile. The occupation distribution 
Psat(dV r S at) is well-approximated by a Poisson distribution 
(e.g., Kravtsov et al. 2004; Yang et al. 2008; Wetzel & 
White 2010), 

\ -^sat p ^ 

P(N sat \M) = (12) 

-i''sat • 

where A = (iV sat |M), so we set {N sat (N sat - 1)\M) = 
(N sat \M) 2 . The two parts of the 1-halo term in equa¬ 
tion m can be thought of as the ‘center-satellite term’ 
and the ‘satellite-satellite term’. 

See Appendix[B]for more details about the modeling of 
the galaxy density profile, halo bias, halo mass function, 
and matter power spectrum. 


3.3. Mock Galaxy Catalogs 

In this section, we describe how we model ‘mock’ 
galaxy catalogs, which is distinct from our analytic model 
and provides an alternative methodology to constrain 
central and satellite galaxy properties and their relations 
with halo mass and halocentric position. The results of 
both models will be presented in Section |4] 

We construct mock galaxy catalogs using the Bolshoi 
dark matter simulation (Klypin et al. 2011 J3 with Rock- 
star phase-space halo-finding algorithm (Behroozi et al. 
2013a, 2013b)EJ The simulation uses a computational 
box with length Lbox = 250 h _1 Mpc and with impressive 
mass resolution (particle mass m = 1.35 x 10 8 h~ 1 M Q ) 
and force resolution (1 /i _1 kpc physical). The model 
for constructing the mocks is described in Skibba et 
al. (2006) and Skibba & Shetlr (2009) with additional 
updates and improvements in the treatment of galaxy 
color distributions and dynamics that will be described 
in Skibba (in prep.). The analytic HOD model (UOMllI) 
assumes the ‘central galaxy paradigm’, in which the most 
massive galaxy is assumed to be the central galaxy of 
a halo and at rest at the halo center, while the mock 


16 Note that Bolshoi has as = 0.82 while the analytic model has 
a s = 0.80. 

11 Note that for consistency with Tinker et al. (2008), we have 
used a M 200 m mass definition for the analytic model, while the 
Rockstar haloes have a Mhnqb mass definition based on the evolv¬ 
ing virial overdensity calibrated by Bryan & Norman (1998). 


catalogs relax these assumption (see Skibba et al. 2011), 
though in practice it has only a minor effect on the small- 
scale projected clustering. Satellite galaxy distributions 
are assumed to follow the same (NFW) density profile as 
in the analytic model, and except for a test below, we 
ignore the subhaloes in the simulation. A version of this 
model was used to construct mock catalogs for Old et al. 
(2015). 

We will assess in Section [4] how consistent the inferred 
HOD parameters as a function of stellar mass from the 
mocks are with the analytic model. In subsequent work, 
we will utilize these models to analyze the clustering 
and distribution of star-forming and quiescent galaxies 
in more detail. 

As an example of our models, we compare a correlation 
function prediction of the HOD-based analytic model and 
mock galaxy catalog in Figure [3] The two models are 
sufficiently consistent where the 1-halo term dominates 
(r p < 500 /i -1 kpc) and where the 2-halo term dominates 
(r p > 3 h _1 Mpc). However, the analytic model’s clus¬ 
tering prediction for the scales in between is too low. 
Observed clustering measurements do exhibit a bump 
in this region (Zehavi et al. 2004), though with a more 
power-law-like behavior than in this analytic model (see 
also Watson et al. 2011). Incorporating scale dependent 
bias and halo exclusion, which reduce the discrepancy 
on these scales, are work in progress (see Appendix iBl). 
although there currently is no ideal way to treat these 
effects or to address the issue of halo truncation at a 
particular radius. 

For comparison and as a test, we also show the clus¬ 
tering prediction of a mock catalog constructed with a 
HOD/SHAM hybrid model like one we developed for 
the comparison project in Knebe et al. (2O15)0- In 
this model, satellite galaxies are directly associated with 
rank-ordered subhaloes and are given their positions and 
velocities, rather than assuming NFW distributions. The 
result is nearly identical to that of our fiducial HOD- 
based mock catalog. At lower stellar masses than we 
probe in this paper, in the regime where a significant 
fraction of satellites are ‘orphans’ whose subhaloes have 
been stripped away (e.g., Wang et al. 2006), galaxies in 
these two types of models have slightly different spatial 
distributions, however (Pujol, Skibba, & Gaztanaga et 
al., in prep.). 

In the next section, we will present results for both sets 
of models, though the constraints on the satellite HOD 
parameters may be more robust in the mocks because 
of their more realistic treatment of the 1-halo to 2-halo 
transition region. 

4. RESULTS: HALO OCCUPATION PARAMETERS 
4.1. Parameter-Fitting Procedure 

For both sets of models, we perforin parameter scans 
and calculate y 2 values, and our procedure yields best-fit 
halo-model parameters and approximate l-a confidence 
intervals. This is not the same as a Monte Carlo Markov 
Chain (MCMC) procedure (see e.g., Tinker et al. 2013), 
though such an analysis should yield nearly identical re- 

18 In practice, it is similar to a full SHAM model in which one 
assumes a stellar-to-halo mass relation and associates central and 
satellite galaxies with haloes and subhaloes in the simulation ac¬ 
cordingly. 














6 


Skibba et al. 



Fig. 3. — Galaxy auto-correlation function comparison for galax¬ 
ies with log(M/M q/i - 1 ) > 11.4 and z ~ 0.7: Bolshoi/Rockstar 
HOD-based mock galaxy catalog (square points), HOD/SHAM 
hybrid mock catalog (diamonds), and HOD-based analytic model 
(black line), including the 1-halo and 2-halo terms (blue and red 
lines). 


suits; our procedure is also similar to that of Wake et al. 

( 2011 ). 

We quantify the total x 2 for a given model and dataset 
with the following: 


N m 


Xtot = J2 *SMF + J2 X* 


i= 1 


(13) 


where Xg MF is for the number density compared to the 
Moustakas et al. (2013) SMF and xt ,* is f° r th e w p {r p ) 
correlation functions compared to M15 at a range of pro¬ 
jected separations (0.1 < r p < 30 /i _1 Mpc). In addition 
to these auto-correlation functions, one could also include 
rank-ordered stellar mass mark correlation functions (see 
Skibba & Sheth 2009; Skibba et al. 2013) when determin¬ 
ing the best-fit models, but such an analysis is beyond 
the scope of this paper. 

In most cases, the host halo masses, n = Mi/M m i n , 
and a are the most robust parameters, while o\ og M and 
Mo are less strongly constrained. (The scatter a\ og M 
is better constrained by other statistics such as satellite 
kinematics and CLFs.) 

The inferred halo-model parameters can change when 
observed abundances (number densities) are incorpo¬ 
rated in the calculation. For the bulk of this paper, we 
present results in which the COSMOS field is included. 
As shown in S14, because of apparently anomalous struc¬ 
tures in this field, it yields clustering strengths and bias 
values much larger than other fields. We split our sam¬ 
ples at the median redshift of z = 0.7 (see Sec. ED, which 
partially splits a large structure in COSMOS (Scoville 
et al. 2013) and reduces this effect. For comparison, in 
Appendix [C] we include the halo-model parameters in¬ 
ferred when the COSMOS field is excluded. In general, 


the results are well within the 1-cr errors of those we ob¬ 
tain here, with the halo masses 0.05-0.10 dex lower in the 
high-redshift range (0.7 < z < 1.2). 

Although we will treat each stellar mass and redshift 
bin independently for the purposes of parameter fitting, 
they are related to each other. For the better constrained 
parameters, fitting formulae as a function of mass and 
redshift can be estimated from the trends we obtain. (See 
Sec. 0 see also Appendix A2 of Skibba & Sheth 2009 for 
luminosity-dependent HOD parameters.) In addition, we 
do not include the covariance matrices of the clustering 
measurements in the model fitting, as they have large 
uncertainties; including the full error covariance matri¬ 
ces, when they are less noisy, versus only the diagonal 
elements usually has only a small effect on clustering 
analyses (Zehavi et al. 2011). We attempt to sufficiently 
finely probe the parameter space so as to estimate 1-cr 
errors of the inferred model parameters. However, there 
are some weak degeneracies between parameters; for ex¬ 
ample, lower Mi/M m j n and higher a both increase the 
satellite fraction (though with different halo mass de¬ 
pendencies) and increase the clustering signal of central- 
satellite and satellite-satellite pairs (though with differ¬ 
ent scale dependencies). 

As a test, we ran models more finely through the pa¬ 
rameter space for the M2 galaxies (M* > 10 10 5 Mq at 
z ~ 0.5), and we show the parameter distributions we 
obtained for M m ; n , Mi/M m i n , and a in Figure|4j. These 
demonstrate the robustness of our results, and the re¬ 
sults shown here are very similar to those we obtain with 
our standard binning in the next section. In particular, 
we obtain logM m i n = 12.132 ± 0.037, /i = 16.80 ± 0.90, 
and a = 1.125 ± 0.053, and in each case the means and 
medians are nearly identical. As stated above, some de¬ 
generacy between the and a parameters can 

be observed. Unless stated otherwise, in the following 
we use halo mass bins of 0.05 dex, Mi/M m ; n bins of 1.0, 
a bins of 0.05, and <J\ og M bins of 0.1. 

4.2. Best-Fitting Models 

Following the procedure in the preceding section, we 
obtain best-fitting halo models with low values of y 2 
relative to the clustering and abundances of galaxies in 
the five catalogs described in Section [2] In Figure [5j 
we show measured projected two-point correlation func¬ 
tions from M15 and best-fitting mock catalogs for galax¬ 
ies with stellar masses log M* > 10.5 at low and high 
redshift {z ~ 0.5 and 0.9). 

In general, we find good agreement, in that our halo 
occupation models reproduce the PRIMUS correlation 
functions from M15 and the cumulative number densities 
(i.e., abundances from the stellar mass functions of Mous¬ 
takas et al. 2013) well. For galaxies with log M* > 10.5 
(middle panels) and with log M* >11.0 (right panels), 
the measured correlation function at z ~ 0.9 is slightly 
higher than at z ~ 0.5, and we infer a slightly higher halo 
mass for the z ~ 0.9 sample. The small-scale correlation 
functions (one-halo terms) are slightly different as well, 
indicating different satellite HOD parameters. 

Next, in Figure |6] we show the mean of the halo occu¬ 
pation distributions, {N ga x) (see Eqns.njj), for the best¬ 
fitting models for the two samples with galaxies of mass 
M* > 10 10 - 5 Mq as an example. Although the shapes 
of the distributions are similar, there are notable dif- 
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log (M mln /h-*M 0 ) Mj/M mln 


Fig. 4. 1-D (diagonal) and 2-D likelihood distributions of the 

model parameters M m ; n , and a resulting from fitting 

the procedure for the M2 galaxies (M* > 10 10 5 Mq at z ~ 0.5). 
Distributions for the top 10% best-fitting models are shown; the 
top 20% have similar distributions. Flat priors were used within 
the ranges shown. 

ferences: in addition to the expected higher halo mass, 
the z ~ 0.9 galaxies also have slightly higher a\ og M and 
slightly lower parameter. 

4.3. HOD Results: inferred model parameters 

We present the mass and redshift-dependent results of 
our analyses in Table [2] and in Figures [7] and [8] The 
host halo masses at a given stellar mass are slightly dif¬ 
ferent when only clustering constraints are used versus 
when those and abundances are used, which is our de¬ 
fault method. For example, for the most massive galaxies 
at * ~ 0.5 (M3), the measured clustering strength implies 
host halo masses of Mh > 10 12 6 M@/h while the number 
density favors a higher halo masses of Mh > 10 12,9 M 0 //i. 
Conversely, for the M* > 10 10 5 Mq galaxies at z ~ 0.9 
(M4), the clustering strength implies host halo masses of 
Mh > 10 12A0 Mq/H while including the number density 
results in a slightly lower value (10 12 25 ), though this may 
be due to anomalously strong clustering in COSMOS. We 
refer the reader to Appendix [Cl for more details. 

We obtain very similar results with the analytic halo 
models (see Table [3]) versus the mock galaxy catalogs. 
The halo masses are nearly identical, while the M\/M m \ n 
values are slightly lower than for the mocks. We have 
tested the issue of the 1-halo to 2-halo transition region 
by artificially inflating the errors of the clustering mea¬ 
surements on those scales (0.5 < r p < 2 /i _1 Mpc), and 
we obtained almost the same HOD parameters (within 
1-cr), although we obtain slightly higher values (by a few 
percent) of the a parameter for the satellite galaxy HOD 
for some of the samples. 

Note that the samples of galaxies selected at or above a 
particular stellar mass have number densities that evolve, 
so that the high-z galaxies are not necessarily progenitors 


of their low-z counterparts (Tojeiro et al. 2012; Leja et al. 
2013). For this reason, some authors have rank ordered 
by luminosity or stellar mass and performed clustering 
analyses at a given number density (e.g., H. Guo et al. 
2013; S14). However, for our range of stellar masses and 
redsliifts, we obtain qualitatively similar results for mass- 
and number density-selected samples. (Note that the 
number densities of M2 and M4 and of M3 and M5 are 
very similar in Table [l]) 

4.4. HOD Results: stellar mass dependence 

First, we show the stellar mass-halo mass relation 
(SHMR) for both redshift ranges (0.2 < z < 0.7 and 
0.7 < z < 1.2) in Figure [7] Here and in what follows, we 
refer to the SHMR as the mean and rms of galaxy stellar 
mass as a function of halo mass, which quantify the dis¬ 
tribution p(M*|M) (e.g., Behrooziet al. 2010; Leauthaud 
et al. 2012). 

Two sets of relations are shown, such that threshold 
stellar masses (M* > 10 10 0 , 10 10 5 , and 10 110 ) are as¬ 
sociated with threshold halo masses, and median stellar 
masses of galaxies in each sample are associated with 
a median halo mass. In both cases, the results are for 
best-fitting mock catalog-based models that minimized 
x 2 (El)- As expected, we find a strong correlation be¬ 
tween the masses at both redshifts, implying a steep slope 
in the SHMR over this mass range. 

Second, we show the satellite galaxy halo occupation 
parameters, Mi/M m ; n and / sa t, in FigurelU The satellite 
fraction is simply defined by the following: 

J Minin dM(dn/dM)(N sat \M) 

/satl *’ f M . dM ( dn/dM) ((N cen \M) + (N sat \M))' 

(14) 

We find that is roughly constant with stellar 

mass, perhaps with some indication of an anticorrelation 
at lower redshift, consistent with H. Guo et al. (2014) 
and Zehavi et al. (2011). The satellite fraction decreases 
rapidly from « 25% at low stellar masses to ~ 10% at 
high masses, consistent with van den Bosch et al. (2007) 
at lower redshift. Most massive galaxies are centrals in 
massive haloes, and they are typically surrounded by less 
massive satellites. (We discuss these issues further in 
Section [5]) 

The values of o\ og M are not strongly constrained, but 
we find that they are approximately constant with mass 
and redshift, usually near a value of ss 0.20. They are 
consistent with the satellite kinematics analysis of More 
et al. (2011), the clustering/lensing analysis of Cacciato 
et al. (2009), and the high-mass galaxy constraints of 
Shankar et al. (2014), while Zheng et al. (2007) and 
Leauthaud et al. (2012) obtained slightly higher values 
of ~ 0.3 and ~ 0.25. Note that the cr\ og M parameter 
contains both a measurement error due to the stellar 
mass measure and redshift error and the intrinsic scatter. 
Moreover, different modeling frameworks (e.g., HOD, 
SHAM, as well as group catalogs and satellite kinemat¬ 
ics) do not necessarily infer the same quantity, highlight¬ 
ing the difficulty of studying it precisely (see Leauthaud 
et al. 2011). 

4.5. HOD Results: redshift evolution 













(h - 1 Mpc) 
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Fig. 5. — Top ten best-fitting halo occupation models (using clustering and abundance constraints) for the mock galaxy catalogs (described 
in Sec. 13.3b and the PRIMUS catalogs (desc ribe d in Sec. [2}. The analytic models produce similar results, though with a 1-halo to 2-halo 
term transition that is too distinct (see Sec. not . Each column shows results for different stellar mass thresholds, and upper (lower) rows 
show results at lower (higher) redshift. Blue points indicate our fiducial measurements from all of the PRIMUS and DEEP2 fields, and 
cyan points indicate our measurements with the COSMOS field excluded. 


TABLE 2 

Halo Occupation Distribution Results: Mock Catalog Analysis 


sample 

logMmin 

(log M h ) 

°"lo gM 


OL 

/sat 

Ml 

11.70±0.05 

12.00±0.10 

0.20±0.10 

17±1 

1.05±0.05 

0.22±0.02 

M2 

12.10±0.10 

12.38±0.10 

0.25±0.10 

17±1 

1.11±0.05 

0.16±0.01 

M3 

12.60±0.10 

12.85±0.10 

0.20±0.10 

16±1 

1.09±0.05 

0.11±0.01 

M4 

12.35±0.10 

12.60±0.10 

0.20±0.05 

14±1 

1.11±0.05 

0.12±0.01 

M5 

12.65±0.10 

12.87±0.10 

0.20±0.10 

15±1 

1.12±0.05 

0.08±0.01 


HOD results for mock galaxy catalogs of the stellar mass and redshift-dependent (volume-limited) PRIMUS catalogs described in 
Table[I] (log Mh) indicates the median halo masses; the mean halo masses are slightly higher (by ~ 0.1 dex in most cases). Ml, M2, 
and M3 are at z ~ 0.5 while M4 and M5 are at z ~ 0.9. Ml has the logM* threshold 10.0, M2 and M4 have 10.5, and M3 and M5 
have 11.0. 


TABLE 3 

Halo Occupation Distribution Results: Analytic Model Analysis 


sample 

logM min 

(log M h ) 

^lo gM 

Mi/M min 

a 

fsat 

Ml 

11.70±0.05 

12.00±0.10 

0.20=t0.10 

18±1 

1.09±0.05 

0.21±0.02 

M2 

12.05±0.10 

12.34±0.10 

0.20=t0.10 

16±1 

1.13±0.05 

0.17±0.01 

M3 

12.60±0.10 

12.85±0.10 

0.20=t0.10 

15±1 

1.10±0.05 

O.lliO.Ol 

M4 

12.35±0.10 

12.60±0.10 

0.20=t0.10 

12±1 

1.10±0.05 

0.13±0.01 

M5 

12.65±0.10 

12.87±0.10 

0.20±0.10 

12±1 

1.03±0.05 

O.lOiO.Ol 


Same as Table l2l but for the analytic halo model of galaxy clustering. 
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Fig. 6. — Mean halo occupation number, ( N\M ) for all, central, 
and satellite galaxies for best-fitting models based on analysis with 
mock catalogs of galaxies with M* > 10 10 , 5 Mq. Low (high) red- 
shift indicated by black (red) lines. The shape of the mean HODs 
of the other three galaxy samples are similar. The best-fit HOD 
parameters of all of the samples are presented in Tables |2] and [3] 



M halo (h-> M 0 ) 


Fig. 7. — Stellar mass-halo mass relation as a function of halo 
mass obtained from the analyses of HOD-based mock catalogs (see 
Table [21. Results are shown for mass thresholds (solid points with 
error bars) and median masses (open points with error bars omit¬ 
ted, for clarity). Low (high) redshift results are indicated by black 
(red) points. 



10 10 10 11 


M* (M 0 ) 

Fig. 8.— A4i/M m ; n parameter (upper panel) and satellite frac¬ 
tion (lower panel) as a function of stellar mass threshold obtained 
from the analyses of HOD-based mock catalogs. Low (high) red- 
shift results are indicated by black (red) points, and stellar mass 
error bars are omitted for clarity. 

We find that the stellar mass-halo mass relation of cen¬ 
tral galaxies appears to evolve with redshift, with a sig¬ 
nificance of 2-3 a. In particular, at stellar masses of 
3 x 1O 1O M0, the SHMR evolves such that a given galaxy 
stellar mass translates to a more massive halo at higher 
redshift. This is approximately consistent with other 
studies in the literature, and we compare and discuss 
them in Section EU 

The AIi/M min parameter also evolves, with a decrease 
at higher redshift of 2.5 and 3.5-cr significance for the 
mock catalogs and analytic model, respectively, and / sa t 
evolves over this range as well. Although haloes with a 
lower value of Mi/M m ; n will host more satellites, note 
that a given stellar mass here does not translate to the 
same halo mass or number density over this redshift 
range, as discussed in the previous section. Haloes at 
fixed mass correspond to rarer density peaks at higher 
redshift and have lower number densities and fewer satel¬ 
lites these earlier epochs. We discuss implications of our 
results for satellite abundances and their formation and 
destruction in Section 15721 The other HOD parameters, 
a and <Ti og M, are approximately constant with redshift. 

In the next section, we interpret these results and dis¬ 
cuss their implications for central and satellite galaxy 
evolution in the context of the halo model. 

5. DISCUSSION: IMPLICATIONS OF THE RESULTS 
5.1. Stellar Mass-Halo Mass Relation 

We now compare our mean stellar-to-halo mass rela¬ 
tion (SHMR) for threshold masses to other studies in 
the literature in Figure El In particular, we compare our 
models of (M|M*) to those of Moster et al. (2010) and 
Behroozi et al. (2010), which involve statistically deter¬ 
mining the SHMR with a given shape and observed abun¬ 
dances (SMFs) as a function of stellar mass and redshift 










10 


Skibba et al. 



Fig. 9. — Stellar mass-halo mass relations obtained in this paper 
at z = 0.5 and 0.9 (as in Fig. [TJ) compared to others in the litera¬ 
ture, where black and red lines indicate lower- and higher-redshift 
results, respectively. Results are shown for Moster et al. (2010; 
dotted lines at z = 0.5 and 0.9); Behroozi et al. (2013; dashed 
lines at z — 0.5 and 0.9), Leauthaud et al. (2012; long-dashed line 
at z = 0.66), Coupon et al. (2015; dot-dashed line at z — 0.83), 
all of which involve fitting relations to binned data, and Zheng et 
al. (2007; triangle points at z ~ 0.95) based on luminosity-binned 
samples. 

and stellar mass-dependent clustering at 2 ~ 0. We also 
compare to Leauthaud et al. (2012), who fit to number 
densities, galaxy-galaxy lensing, and angular clustering 
in COSMOS, and to Coupon et al. (2015), who fit to 
number densities, lensing, and clustering in CFHTLenS. 
Necessary adjustments have been made to their quoted 
results when different IMFs or Hubble constant conven¬ 
tions were assumed. Other relevant results include lumi¬ 
nosity dependent clustering analyses with VIPERS (de 
la Torre et al. 2013) and the COSMOS and CFHTLenS 
results for red and blue galaxies (Tinker et al. 2013; Hud¬ 
son et al. 2015), though as these studies do not constrain 
the exact relation being studied here, we cannot com¬ 
pare our results directly. Within the l-cr error bars, our 
2 ~ 0.5 results are consistent with most of these studies, 
though there are some notable differences: for example, 
we obtain larger halo masses than Leauthaud et al. for 
galaxies with M* ~ 10 1O - 5 M 0 . 

At lower stellar masses, our results suggest a slightly 
stronger redshift evolution than obtained by others. For 
example, at M* = 1O 1O5 M0, we estimate an increase 
in halo mass by 0.25-0.30 ± 0.10 dex from 2 ~ 0.5 to 
2 ~ 0.9, while the models of Moster et al. and Behroozi 
et al., which are not constrained by high-redsliift cluster¬ 
ing, predict a mass increase of no more than 0.1 dex over 
this redshift range—less than half as much evolution as 
we find. Our estimated SHMR evolution is not statis¬ 
tically significant (0.15 ± 0.10 dex) when the COSMOS 
field is excluded (Tabled); however, the analysis of Leau¬ 
thaud et al. finds less evolution within COSMOS based 
on the mean stellar and halo masses. In another lensing 
analysis, Hudson et al. (2015) find no SHMR evolution 


for blue galaxies and ss 0.2 dex evolution for red ones. 

Because of our limited dynamic range, it is difficult 
to fit a function to our SHMR results. Common double 
power-law parameterizations include the following ones 
of Moster et al. (2010) and Behroozi et al. (2010), re¬ 
spectively: 



log(H4(M*)) =log(Mi) + /?log(——) 

- M *,0 

(M*/M*, „)* 1 

+ 1 + (M,/Af, i0 )-T 2’ 


(16) 


where /? and 5 quantify the low- and high-mass slopes. 
Moster et al. obtain low- and high-mass slopes of ap¬ 
proximately 2.2 and 0.5, respectively, while Behroozi et 
al. obtain slopes of approximately 2.6 and 0.3. It is pos¬ 
sible to measure the relation’s approximate slope to our 
model results as well, though the precise transition from 
the low- to high-mass regime is not clear. For galaxies at 
2 ~ 0.5 and using the threshold masses, if M* oc Mj®, 
then /3 ~ 1.3 at Mh ~ l0 12 M@/h. and ft ss 1.0 at 
Mh ~ 10 12 6 Mq//i. When the median masses are used, 
the slope is slightly shallower by ~ 0.1. At 2 ~ 0.9, 
the slope is slightly steeper than that of other models 
for halo masses less than Mh ~ 10 12 ' 6 Mq//i. At higher 
masses, for group and cluster haloes, the Moster et al. 
and Behroozi et al. models predict that the slope in¬ 
creases by about 20-25% between 2 ~ 0.5 and 0.9, though 
our model constraints are not sufficiently precise at such 
masses to determine whether we obtain similar evolution. 

Rather than analyze this relation, some authors take 
the ratio, /* = which does not contain new in¬ 

formation but more clearly indicates the peak or ‘pivot’ 
mass scale (e.g., van den Bosch et al. 2007; Leauthaud et 
al. 2012). Based on our results in Tableland Figure [3 
we find that /* peaks at a mass scale of 12.0 < log Mh < 
12.4 at a value of 0.028, consistent with the literature 
and well below the universal baryon fraction of approxi¬ 
mately 15.5% (assuming VL^h 2 = 0.022 and fl c h 2 = 0.120 
from Planck collaboration et al. 2014; see also Fukugita 
& Peebles 2004). At 2 ~ 0.9, the pivot (log) halo mass 
appears to shift upward to « 12.3-12.6, slightly higher 
than studies in the literature, though a wider dynamic 
range with high completeness would be necessary to an¬ 
alyze this further and obtain more precise results. 

It is well known that the SHMR has a significant 
amount of scatter (More et al. 2011; Skibba et al. 2011; 
Rodrfguez-Puebla et al. 2013). The relation’s scatter in 
our best-fit models, quantified by cr\ og M, is consistent 
with the literature but has large error bars. Note that 
our central galaxy HOD assumes a lognormal distribu¬ 
tion of stellar mass at fixed halo mass (Eqn.[3|; although 
this assumption is consistent with the current data avail¬ 
able, it may not be entirely accurate. Also note that 
SHMR studies in the literature do not all assume the 
some halo mass definition (see App. El, though even for 
the same definition, halo occupation and semi-analytic 
models make a wide range of predictions for the SHMR 
at a given halo mass (Knebe et al. 2015). 
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Finally, mean or median SHMRs of central and satel¬ 
lite galaxies are similar when each object’s halo mass (or 
circular velocity) at the approximate time of accretion is 
used, though we refer the reader to Rodriguez-Puebla et 
al. (2012, 2013) and Watson & Conroy (2013) for anal¬ 
yses of some differences between them. However, as a 
function of parent halo mass, satellites have surprising 
mass distributions and relations, with implications for 
their co-evolution with the subhalos that host them, and 
we turn to these issues in the following subsection. 

5.2. Satellite and Subhalo Abundances 

In Section 14.51 we noted that the satellite HOD (de¬ 
scribed in Sec. m , and in particular the pa¬ 

rameter (or /i) appears to evolve, consistent with other 
studies (e.g., de la Torre et al. 2013). Combining our 
satellite HOD results with lower-redsliift results from 
the SDSS (Zehavi et al. 2011) and higher-redshift re¬ 
sults at 2 ~ 1 (Zheng et al. 2007) and at z ~ 1.5 (Wake 
et al. 2011; Martinez-Manso et al. 2014) imply a clear 
redshift dependence. For example, applying the func¬ 
tion Mi/M min (z) = ho — (3z to all of these results yields 
Ho ~ 19-20 and /3 « 6-7. 


We can use our results to determine the approximate 
amount of stellar mass in a diffuse stellar halo or the intr¬ 
acluster medium, likely due to disrupted satellite galax¬ 
ies. Such a calculation is made by comparing the halo 
mass fraction in subhalos to the stellar mass fraction in 
satellites (Skibba et al. 2007; White et al. 2007; Yang 
et al. 2009). We utilize our constraints on the stellar 
mass and redshift dependent halo occupation parame¬ 
ters, especially and we compare the relative 

abundances of satellite galaxies to that of subhaloes esti¬ 
mated from the subhalo mass function (MF). From Gio- 
coli et al. (2010), the subhalo MF can be expressed as 
the following: 


dN(m\M) 

dm 


dm = N(z) 


M 


0.9 


exp 


10 12 h~ 1 M l 

-b( 


m 

\M 


o.i 

dm 

m 


(17) 


where m is the subhalo mass at the redshift of interest, 
N(z) « 0.0148(1 + z) 1 / 2 and /? ~ 12.2715. If we neglect 
the exponential in (Eon. [171) . then 


5.2.1. Satellite Mass Function 

The decreasing with increasing redshift is 

related to the competition between merging and satellite 
destruction (see Conroy et al. 2006; Wetzel et al. 2009). 
These results have implications for satellite abundances 
and fractions as a function of mass and redshift. One 
way to study this is to revisit the analysis of satellite 
galaxy and subhalo occupation distributions of Skibba 
et al. (2007); we refer the reader to that paper for more 
details. 

Firstly, note that although the most massive satellite’s 
mass or luminosity will scale with the host halo, the mean 
satellite mass or luminosity is nearly flat over many or¬ 
ders of magnitude of parent halo mass (see also Skibba 
et al. 2011; Paranjape & Slieth 2012). Our analysis also 
predicts that not just the mean, but the shape of the 
satellite galaxy stellar mass function is approximately 
independent of halo mass. 

As we argued in Skibba et al. (2007) for galaxy lumi¬ 
nosities, the power-law shape of the mean satellite HOD 
|3|) and the approximate exponential shape of the stellar- 
to-halo mass relation (except at low masseJ^ll implies 
that the satellite conditional SMF, p(M* iSat |M), is ap¬ 
proximately has a Schechter-likc function. This is borne 
out in the analysis of central and satellite conditional 
stellar mass distributions in group catalogs by Yang et 
al. (2009). The evolving h parameter with redshift only 
implies that the amplitude of the satellite conditional 
SMF evolves, while the shape of it remains constant at 
least out to z ~ 0.9. Therefore, satellite abundances or 
occupation numbers of haloes evolve, but low and high- 
redsliift haloes have similar relative abundances of low 
vs high-mass satellites. 


N{> m\M ) = 



dN(m\M) 

dm 


dm 


N(z) / M \ 
0.9 V10 



If we use Mi to denote the value of M at which the 
number of subhaloes is unity, then the expression above 
implies that 



60.8 /10 12 /i -1 Mq 

(1 + z ) 1 / 2 \ m 


(19) 


Going a step further, using the subhalo MF in fEan.fTTl 
and neglecting the exponential), the mass fraction in 
subhaloes is approximately given by 


r M 


/sub(M) = 


m dN[m\M) 


M 


dm 


dm 


= 0.148(1 + z) l/2 ^ 


M 


10 


© 


0.1 

> ( 20 ) 


where if stars only form in sufficiently massive objects, 
the lower limit to this integral may be greater than zero. 
The p = parameter of galaxies is as low as 

« 12 at z ~ 1, consistent with Martinez-Manso et al. 
(2015) and other studies, but for subhaloes, the ratio is 
much larger, implying that some mass has been stripped 
from satellites, presumably contributing to stellar mass 
in the intracluster light (ICL). Revisiting the calculation 
in Skibba et al. (2007) but with new clustering and sub¬ 
halo constraints, we argue that the mass fraction in the 


5.2.2. Intracluster Mass 


19 Over a wide mass range, including masses well below the knee 
of the SMF, a double power-law shape is more accurate (Moster et 
al. 2010; Behroozi et al. 2010; see also Moustakas et al. 2013). 


20 According to Giocoli et al. (2010), neglecting the exponen¬ 
tial is not entirely accurate, but it is sufficiently accurate for the 
very approximate calculation here. The integral can be performed 
analytically with the exponential and yields incomplete gamma 
functions in the result. 
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intracluster medium is given by 

/ IC m(M,M.) = 1 - (/subtM.z) 60 ' 8 ^!^ V2 ) ( 21 ) 

<AU|M,M*) + (M*, cen |M)/(M*, sat |M)\ 

(AT sat |M,M*) /’ 

where the halo redshift dependencies cancel, imply¬ 
ing that this fraction’s evolution depends primarily on 
—that is, on how the halo occupation number of 
subhaloes relative to satellite galaxies evolves. 

The halo model predicts the ICL fraction to increase 
with host halo mass, in agreement with other analyses 
(Murante et al. 2007; Purcell et al. 2007), which also 
highlight uncertainties about the fate of ‘orphan’ satel¬ 
lites over time (Conroy et al. 2007; Contini et al. 2014). 

For halo masses of M ~ 10 14 /i _ 1 M 0 and galaxies with 
M* > lO la 5 M 0 , we estimate that this fraction is ap¬ 
proximately 5% at z ~ 0.5 and less than 2% at z ~ 0.9. 

For halo masses of M ~ 10 15 /i - 1 M@ the ICL fraction 
is approximately three times larger, though the rarity of 
such massive haloes at high redshift makes this calcula¬ 
tion even more uncertain. It appears that the stripping 
and disruption of satellites indeed contribute to the ICL 
over the 2.5 Gyr between these epochs, but additional 
analysis with more precise halo-model parameters is re¬ 
quired to study this further. 

6. CONCLUSIONS 

In this paper we analyze the stellar mass and redshift 
dependent spatial clustering and abundances of galaxies 
in PRIMUS and DEEP2 at 0.2 < z < 1.2 with ACDM 
halo models of galaxy clustering. In order to obtain ro¬ 
bust results for model parameters, we perform the anal¬ 
ysis with two independent sets of models: HOD-based 
mock galaxy catalogs and analytic models. 

We summarize our main conclusions as follows: 

• Both sets of models are able to accurately repro¬ 
duce the stellar mass and redshift dependent corre¬ 
lation functions and number densities. The best-fit 
models yield consistent and robust halo occupation 
parameters for central and satellite galaxies. 

• The results for central galaxies constrain the 
redshift evolution of the mean stellar-to-halo 
mass relation (SHMR). For halo masses below 
10 126 /i _ 1 Mq, we find that the SHMR appears 
to evolve significantly, such that galaxies of a 
given stellar mass are associated with more mas¬ 
sive haloes at higher redsliifts. We find that the 
ratio M*/Mh a i 0 peaks at a value of 0.028 at the 
mass scale 12.0 < log(Mh a i 0 //i _ 1 MQ) < 12.4. We 
do not obtain statistically significant constraints on 
the evolution of the scatter in the SHMR. 

• The satellite fraction increases rapidly with de¬ 
creasing stellar mass and the Mi/M m i n parame¬ 
ter, which quantifies the critical mass above which 
haloes host at least one satellite, is more redshift 
than mass dependent, with lower values (~ 12-15) 
at high redshift. We use the HOD of satellites 
and subhaloes to estimate rough constraints on the 
mass fraction of disrupted satellites that contribute 
to the intracluster medium. 


Our joint analysis with analytic models and mock cat¬ 
alogs, which produce consistent results, demonstrates 
their robustness and the strength of our conclusions. 
Abundance-matching with numerical simulations and 
studies with hydrodynamic simulations are currently 
popular in the field, but analytic models remain an im¬ 
portant tool as well. All models have some assumptions, 
uncertainties, and free parameters, and every model has 
advantages and shortcomings; therefore, it is important 
to utilize multiple types of models when possible. The 
advantages of analytic models include their speed and 
computational inexpensiveness and their ability to pa¬ 
rameterize physical processes and galaxy-halo relations 
in a way that aids understanding their origins. 

By comparing galaxy counterparts at different red- 
shifts, one may constrain the extent to which galaxies of 
a given mass or number density grow by in situ star for¬ 
mation vis-a-vis mergers with neighboring galaxies over a 
given redshift range (Wake et al. 2008; Zehavi et al. 2012; 
Lackner et al. 2012). For example, Wake et al. (2008) uti¬ 
lized a similar analytic halo model at 0.19 < 2 < 0.55 as 
the one used here to estimate the central galaxy merger 
rate of luminous red galaxies. For the galaxy samples 
used here, however, the clustering measurements and 
abundances are consistent with passive evolution, sug¬ 
gesting a very low merger rate, though the uncertainties 
are too large to make quantitative estimates (S14). Nev¬ 
ertheless, Seo et al. (2009) argue that a shoulder in the 
HOD (Mi M m i n ) implies that a galaxy population has 
not undergone passive evolution, and this remains an im¬ 
portant issue for future research. In addition, analysis of 
the clustering of star-forming and quiescent galaxies out 
to z ~ 1 in PRIMUS and DEEP2 would greatly benefit 
the field and would constrain the quenched fractions of 
central and satellite galaxies as a function of mass and 
redshift, which would be important for distinguishing be¬ 
tween competing models (Cohn & White 2013). 

Furthermore, as described in Section [3l throughout 
this paper we have assumed that the evolving spatial 
distributions of galaxies and their correlations with the 
dark matter haloes are primarily determined by the mass 
of the haloes. This implies an implicit assumption that 
‘galaxy assembly bias,’ where galaxies’ distributions de¬ 
pend on the assembly history of systems at fixed halo 
mass, is such a small effect as to be negligible for the 
analysis. Some recent studies in the literature (Zentner 
et al. 2014; Hearin et al. 2014) argue that this assump¬ 
tion may be incorrect and therefore inferred model pa¬ 
rameters and galaxy-halo correlations, such as the ones 
obtained in this paper, may be biased. Addressing this 
question is beyond the scope of this paper, but it too will 
be the focus of subsequent research. 

Finally, this paper is complementary to other current 
and upcoming work. In particular, stellar mass and SFR 
dependent clustering measurements in PRIMUS will be 
presented in Mendez et al. (in prep.; M15), luminosity 
and color dependent cross-correlation functions are pre¬ 
sented in Bray et al. (2015), and the environmental de¬ 
pendence of stellar mass functions are in Hahn et al. 
(2015). In addition, complementary modeling of SFR 
dependent clustering will be presented in Watson et al. 
(in prep.) extending the work in Watson et al. (2015), 
and details about the models used for constructing mock 
galaxy and group catalogs will be in Skibba (in prep.). 
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APPENDIX 

HALO MASSES AND HALO FINDING ALGORITHMS 

As noted in Section [3j throughout this paper we use halo masses and radii defined using a virial overdensity 200 
times the mean density of the Universe for the analytic model or defined using the evolving virial overdensity specified 
by Bryan & Norman (1998) for the mock catalogs. These and other halo definitions typically used in the literature, 
including spherical overdensities (SO), Friends-of-Friends (FoF), and phase-space ones, will have only minor effects on 
the kinds of quantitative results presented here, and the qualitative trends will remain unchanged. 

In general, SO halo finders tend to impose a more spherical geometry on the resulting systems, while FOF sometimes 
links neighboring objects via tenuous bridges of particles. In either case, the choice of virial overdensity and halo 
membership can be important. Dynamically unrelaxed haloes, which are common, as well as poorly resolved haloes 
have poorly estimated or biased masses, concentrations, and other affected parameters (e.g., Skibba & Maccio 2011). 
In addition, some authors adopt different circular velocity definitions, such as Umax and Vp e ak, which may affect some 
model-dependent results (Reddick et al. 2013; Behroozi et al. 2014). 

A detailed analysis of halo definitions is beyond the scope of this paper but will be important for advancing the 
fields of galaxy formation and large-scale structure formation. For more studies of the effects of halo definitions, we 
refer the reader to Knebe et al. (2011), Zemp (2014), Klypin et al. (2014); and More, Diemer & Kravtsov (2015); for 
studies of subhalo definitions, we refer the reader to Onions et al. (2012) and Pujol et al. (2014). 

HALO MASS FUNCTION AND BIAS 
Matter Power Spectrum 

The linear matter power spectrum is described by 

P" n m (k,z) ex T\k)k n °D\z)/D\ 0) (Bl) 

where we use the Smith et al. (2003) matter power spectrum model and the Eisenstein & Hu (1998) transfer function, 
which assumes a CMB temperature of 2.725 K. D(z ) is the linear growth factor and we have set the spectral index 
n, = 1 (see e.g., Dodelson 2003). 

The mass variance is described by 
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which is set such that it is equal to erf for R = 8h *Mpc at z = 0. IE is the Fourier transform of the top-hat window 
function of radius R: 

W(kR) = 3|sin(tfl) - ° 8(tfl) 1 (B3) 

(.kRy 


Halo Mass Function 

We assume the Tinker et al. (2008) halo mass functiorp 1 ! in which haloes are identified with a SO algorithm with 
A = 200 and p m (z) = VL m {z) p cr it{z) = p m (0)(l + z) 3 . Halo abundances are assumed to follow a universal function in 
terms of the mass fraction of matter in peaks of height v = 5 c /a(M, z ), where 8 C = 0.15(127r) 2//3 ~ 1.68(0. The halo 
abundances are well described by the following functional form: 


where 
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dM 
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(B4) 

(B5) 


and the parameters A, a, and b are redsliift dependent and constrained by numerical simulations that are described 
in Tinker et al. (2008). 

We find that this analytic function yields halo abundances as a function of mass and redshift that are approximately 
consistent with Bolshoi (Klypin et al. 2011) Rockstar (Behroozi et al. 2013) catalogs in which masses are defined with 
an evolving virial overdensity (Bryan & Norman 1998). Note that the halo finder used by Tinker et al. allows for 
overlapping haloes. 

The characteristic mass of the halo mass function is defined at the scale at which a typical peak {v = 1) collapses 
at a given redshift: a(M*,z) = a(M* 1 0)D(z) = 6 c (z). M*(z) decreases with increasing redshift, from a few times 
10 12 h~ 1 M Q at z ~ 0 to below 10 7 h~ 1 M Q at 2 ~ 6 (Mo & White 2002). 


Halo Bias Function 

We adopt the Tinker et al. (2010) halo bias function with A = 200, in which the bias is expressed in the following 
flexible form 

b(v) = 1 - A V f + Bv b + Cu c , (B6) 

w v a + <5“ v ’ 

where the parameters A, a, and C are a function of A. This bias function is forced to obey the relation 


J dv b{v)f(v) = 1, 


(B7) 


where f(v) is the halo mass function [such that f(a) = vf(v) using the notation in (IB5I) ]. The Tinker et al. (2010) 
departs from the bias model of Sheth, Mo, & Tormen (200l) at very low and very high values of v. We find that the 
alter nate bias function of Tinker et al. (2010) using the peak-background split yields similar clustering predictions as 
dBTill . 

Scale-dependent bias (at fixed mass) and halo exclusion effects could slightly affect the transition from the 1-halo to 
2-halo terms and the inferred satellite occupation parameters, and an analysis of models of these effects is the subject 
of ongoing work, one approach is that of van den Bosch et al. (2013) and Tinker et al. (2005): 


where 
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(B9) 


Anin (Al) is the cutoff that accounts for halo exclusion, and £(r) quantifies the scale dependence of bias based on some 
fitting functions. 

Considering that the matter power spectrum is modeled in Fourier space, we have opted to model galaxy clustering 
in Fourier space as well. Therefore, implementing the formulation above for our analytic models would require an 
additional integral and would be more computationally expensive. Another promising approach of scale-dependent 


21 This is a modification from our previous models, in which we 
used the Sheth, Mo, & Tormen (2001) mass function. 


22 Note that van den Bosch et al. (2012) use a slightly different 
definition that includes a factor of [Q m (z)] 0 ' 0055 /Z)( 2 ). 
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bias is the excursion set analysis of Musso et al. (2012) and Paranjape et al. (2013) (cf., Smith et al. 2007; Desjacques 
et al. 2010). This approach predicts the following form for the bias factors: 


b 


n 



(B10) 


where we are interested only in the n = 1 case, and 


S X =J dink A 2 (k)W(kR)W(kRo) 

e x =2 d\nS x /dins (Bll) 


are cross-correlations between the mass overdensity field smoothed on the large scale Rq and the Lagrangian scale of 
the halo R. 

More work is required to accurately model scale-dependent bias and halo exclusion. We also find degeneracies 
between the parameters of these models and between them and the treatment of halo profile truncation (see the 
following section). 


Density Profile 


For the halo density profiles, we assume spherical symmetry and adopt a Navarro, Frenk & White (1997; NFW) 
profile: 


P(r\ M ) 


_ Ps _ 

(r/r s )(l + r/r s ) 2 


(B12) 


an NFW profile is not a sufficiently accurate description of the small-scale distribution of satellite galaxies (Watson et 
al. 2012; Piscionere et al. 2014), but this affects the clustering at smaller separations than we can accurately probe. 
The Fourier transform of the NFW density profile is the following: 


u(k\M) 



dr 477 r 2 


sin (kr) p(r\M) 
kr M 


(B13) 


For an NFW profile, the integral can be computed analytically (see also Scoccimarro et al. 2001; Cooray & Sheth 
2002), and the solution to the indefinite integral is 


u(k\M) 


/(c) 


sin(fcr s )Si[fc(r + r s )] + cos(fcr s )Ci[/c(r + r s )] 


sin(fcr) 
r + r s 


(B14) 


where Ci(cc) and Si(:r) are the cosine and sine integrals, respectively, c = r v \ T /r s and the virial radius is defined to be 
r 2 oo m- The expression in (IB14I) is evaluated from r = 0 out to r = r v j r , but if the profile is not truncated at the virial 
radius, then the above reduces to /(c)sin(Ar s )7r/2 as r —> oo (see Sheth et al. 2001 for a discus sion of halo truncation). 
Other profiles such as a Moore et al. (1999) or Einasto (1965) profile may be substituted in (IB 131) . but they must be 
integrated numerically. 

We assume the redshift-dependent concentration-mass relation of Munoz-Cuartas et al. (2011), and in our mock 
catalogs, we include the scatter in the relation. The relation is fitted for dynamically relaxed haloes identified with a 
spherical overdensity algorithm and has the form 

log(c) = a(z)log(M vir /[/i -1 M Q ]) + b(z) (B15) 

We obtain similar clustering predictions using the concentration mass relation of Maccio et al. (2008). 

The concentration and velocity dispersion of galaxies and subhaloes is different than that of dark matter particles 
(e.g., Munari et al. 2013; Old et al. 2013), and we attempt to account for this effect as well. Note that studies in the 
literature currently disagree about the extent to which satellite galaxies are less concentrated than dark matter (see 
also Yang et al. 2005; Hansen et al. 2005). 


COSMIC VARIANCE AND THE EFFECT OF THE COSMOS FIELD 

In Section 7.5 of S14, we discussed the effects of ‘cosmic variance,’ or field-to-field fluctuations, and how PRIMUS’s 
seven science fields help to assess them. For high-redshift galaxy samples, we obtained statistically significant differences 
(up to 3 -a) between the clustering signal within COSMOS and other fields, though the samples were selected differently 
in S14 than in this paper. 

In Table HI we list inferred HOD parameters for the “nocosmos” samples in M15, which are like the default ones 
but exclude the COSMOS field, which has multiple rare large structures. Compared to the results in Section U the 
inferred halo-model parameters, especially at high redshift, are slightly different here than in Tabled For example, 
the inferred high-redshift halo masses are approximately 20% (0.1 dex) lower when COSMOS is excluded. The other 
model parameters are very similar. 
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TABLE 4 

Halo Occupation Distribution Results Excluding COSMOS Field 


sample 

logM min 

(log M h ) 

^logM 

M 1 /M min 

a 

/sat 

Ml 

11.70±0.05 

12.00±0.10 

0.20=t0.10 

17±1 

1.09±0.05 

0.22±0.02 

M2 

12.10±0.10 

12.38±0.10 

0.20=t0.10 

17±1 

1.11±0.05 

0.16±0.01 

M3 

12.60±0.10 

12.85±0.10 

0.20=t0.10 

16±1 

1.10±0.05 

O.lliO.Ol 

M4 

12.25±0.10 

12.53±0.10 

0.20±0.10 

13±1 

1.10±0.05 

0.13±0.01 

M5 

12.60±0.10 

12.83±0.10 

0.20±0.10 

15±1 

1.07±0.05 

0.09±0.01 


HOD results with mock galaxy catalogs for PRIMUS catalogs with the COSMOS field excluded. These results are analogous to 
those presented in Table 121 

But what does it mean to exclude or include a region with seemingly anomalous structures (or voids)? For a 
discussion of this issue, see Norberg et al. (2011) and Zeliavi et al. (2011). One could attempt to assess the effects of 
COSMOS’S structures while using mocks and numerical simulations that include a regions with various large structures, 
but that would require better knowledge of the mass and spatial extent of the structures in the field and is beyond the 
scope of this paper. 






